

In [1]:

    
%pylab inline
from scipy.integrate import simps









    



Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.



In [16]:

    
%run ../paramless.py



In [3]:

    
from JSAnimation.IPython_display import display_animation

Evolution of metabolic investment strategies

Here I test the function simulation framework with the example of section 6 in Dieckmann et al. JTB (2006)

We will evolve a function $x(a)$, where $a \in [0,1]$ (s.t $x(a) > 0 \forall a$ )

Define fitness and mutations

The domain of $x(\cdot)$ is $(0,1)$.



In [4]:

    
domain = np.arange(0.001, 1.0, 0.01)

The energy is defined as $E(x) = \int x(a) da$, which for simulation purposes we will compute using simpson's.



In [5]:

    
def energy(surface, domain, **kwargs):
    return simps(surface, domain)

Metabolic efficiency is defined as $e(a, x(a)) = \frac{x(a)}{x(a)+a}$. Then, the total intake of resources is given by $G(x) = \int r(a)e(a, (a)) da$



In [6]:

    
def resource_density(a, **kwargs):
    return 4.0 * (1.0 - a)

def total_intake(surface, domain, **kwargs):
    #create a list of pairs x(a), a for a in domain
    tuples_surface_domain = zip(surface, domain)
    #r(a)*(x(a)/x(a)+a) for all a in domain
    integrand = [resource_density(surface_domain[1])*(surface_domain[0]/(surface_domain[0]+surface_domain[1])) for surface_domain in tuples_surface_domain]
    #integrate over domain
    return simps(integrand, domain)

With the total intake and the energy we are ready to compute fitness.



In [7]:

    
def metabolic_fitness(resident, mutant, c, domain, **kwargs):
    fitness_resident = total_intake(resident, domain, **kwargs) - c*energy(resident, domain, **kwargs)
    fitness_mutant = total_intake(mutant, domain,**kwargs) - c*energy(mutant, domain, **kwargs)
    return fitness_resident, fitness_mutant

Let's evolve things

We use the point mutation that does not allow for negative values and conserves energy.



In [8]:

    
number_of_generations=5000
initial_surface = np.zeros_like(domain)

First with soft mutations



In [9]:

    
ans, series = evolve(initial_surface, 
                     fitness_function=metabolic_fitness, 
                     mutation_function=gaussian_mutation,
                     iterations=number_of_generations,
                     domain=domain,
                     mutation_epsilon=0.01,
                     c=0.5, 
                     return_time_series=True,
                     atol=1e-8, width=0.02, lower_bound=0.0)



In [10]:

    
plot(domain, ans)









    Out[10]:





[<matplotlib.lines.Line2D at 0x107f0ced0>]



In [11]:

    
def prediction_investment(a):
    c=0.5
    if resource_density(a) > 0.5*a:
        return math.sqrt((a/c)*resource_density(a))-a
    return 0.0



In [12]:

    
target = [prediction_investment(a) for a in domain]



In [13]:

    
plot(domain, target)
plot(domain, ans)









    Out[13]:





[<matplotlib.lines.Line2D at 0x107e946d0>]



In [14]:

    
ani = create_video_from_time_series(series, target, domain, filename='/Users/garcia/Desktop/metabolic_soft.mp4', approximate_number_of_frames=300, record_every=50)
display_animation(ani)

Now with point mutation



In [17]:

    
number_of_generations=50000
ans_stiff, series_stiff = evolve(initial_surface, 
                     fitness_function=metabolic_fitness, 
                     mutation_function=point_mutation,
                     iterations=number_of_generations,
                     domain=domain,
                     mutation_epsilon=0.01,
                     c=0.5, 
                     return_time_series=True,
                     atol=1e-8, lower_bound=0.0)



In [18]:

    
plot(domain, target)
plot(domain, ans_stiff)









    Out[18]:





[<matplotlib.lines.Line2D at 0x1086105d0>]



In [19]:

    
ani_stiff = create_video_from_time_series(series_stiff, target, domain, filename='/Users/garcia/Desktop/metabolic_soft.mp4', approximate_number_of_frames=100, record_every=1000)
display_animation(ani_stiff)